Pathway Heatmaps (Enrichr)
# Get enrichr results and save to file
filename <- file.path(here::here(), "analysis", "avg-expr-GSE137001", "enriched_res.Rdata")
padj_cutoff <- 0.01
lfc_cutoff <- 0.5
dbs <- c("GO_Biological_Process_2021")
if (!file.exists(filename)) {
enriched_res <- lapply(lfc_dfs, function(df) {
genes <- df %>%
filter(padj <= padj_cutoff & abs(log2FoldChange) >= lfc_cutoff) %>%
pull(SYMBOL)
return(enrichR::enrichr(genes, dbs))
})
save(enriched_res, file = filename)
} else {
load(filename)
}
Combined Score
# Extract combined scores from each enrichr result
pathways <- c("GO:0009267", "GO:0034198", "GO:0036003", "GO:0036500", "GO:1990928", "GO:0034599", "GO:0034976", "GO:0016236", "GO:0010506")
combined_scores <- lapply(enriched_res, function(res) {
c_score <- res[["GO_Biological_Process_2021"]] %>%
filter(str_detect(Term, paste(pathways, collapse="|"))) %>%
arrange(Term) %>%
pull(Combined.Score)
return(c_score)
})
combined_scores <- as.data.frame(combined_scores)
row.names(combined_scores) <- enriched_res[[1]][["GO_Biological_Process_2021"]] %>%
filter(str_detect(Term, paste(pathways, collapse="|"))) %>%
arrange(Term) %>%
mutate(Term = str_wrap(Term, width = 60)) %>%
pull(Term)
pheatmap::pheatmap(
combined_scores,
cluster_cols = FALSE,
scale = "row",
)

-log10(adjusted p-value)
# Extract combined scores from each enrichr result
pathways <- c("GO:0009267", "GO:0034198", "GO:0036003", "GO:0036500", "GO:1990928", "GO:0034599", "GO:0034976", "GO:0016236", "GO:0010506")
combined_scores <- lapply(enriched_res, function(res) {
c_score <- res[["GO_Biological_Process_2021"]] %>%
filter(str_detect(Term, paste(pathways, collapse="|"))) %>%
arrange(Term) %>%
mutate(neglog10 = -log10(Adjusted.P.value)) %>%
pull(neglog10)
return(c_score)
})
combined_scores <- as.data.frame(combined_scores)
row.names(combined_scores) <- enriched_res[[1]][["GO_Biological_Process_2021"]] %>%
filter(str_detect(Term, paste(pathways, collapse="|"))) %>%
arrange(Term) %>%
mutate(Term = str_wrap(Term, width = 60)) %>%
pull(Term)
pheatmap::pheatmap(
combined_scores,
cluster_cols = FALSE,
scale = "row",
)

DEGs and Volcano
state_d2_vs_MEF
res <- lfc_dfs[[x]]
padj_cutoff <- 0.01
lfc_cutoff <- 0.5
res %>%
select(SYMBOL, log2FoldChange, padj) %>%
filter(padj <= padj_cutoff & abs(log2FoldChange) >= lfc_cutoff) %>%
arrange(padj) %>%
DT::datatable() %>%
DT::formatSignif(columns = 2:3, digits = 3)
EnhancedVolcano::EnhancedVolcano(
res,
lab = res$SYMBOL,
x = "log2FoldChange",
y = "padj",
pCutoff = padj_cutoff,
FCcutoff = lfc_cutoff,
)

state_d4_vs_MEF
res <- lfc_dfs[[x]]
padj_cutoff <- 0.01
lfc_cutoff <- 0.5
res %>%
select(SYMBOL, log2FoldChange, padj) %>%
filter(padj <= padj_cutoff & abs(log2FoldChange) >= lfc_cutoff) %>%
arrange(padj) %>%
DT::datatable() %>%
DT::formatSignif(columns = 2:3, digits = 3)
EnhancedVolcano::EnhancedVolcano(
res,
lab = res$SYMBOL,
x = "log2FoldChange",
y = "padj",
pCutoff = padj_cutoff,
FCcutoff = lfc_cutoff,
)

state_d6_vs_MEF
res <- lfc_dfs[[x]]
padj_cutoff <- 0.01
lfc_cutoff <- 0.5
res %>%
select(SYMBOL, log2FoldChange, padj) %>%
filter(padj <= padj_cutoff & abs(log2FoldChange) >= lfc_cutoff) %>%
arrange(padj) %>%
DT::datatable() %>%
DT::formatSignif(columns = 2:3, digits = 3)
EnhancedVolcano::EnhancedVolcano(
res,
lab = res$SYMBOL,
x = "log2FoldChange",
y = "padj",
pCutoff = padj_cutoff,
FCcutoff = lfc_cutoff,
)

state_d9_vs_MEF
res <- lfc_dfs[[x]]
padj_cutoff <- 0.01
lfc_cutoff <- 0.5
res %>%
select(SYMBOL, log2FoldChange, padj) %>%
filter(padj <= padj_cutoff & abs(log2FoldChange) >= lfc_cutoff) %>%
arrange(padj) %>%
DT::datatable() %>%
DT::formatSignif(columns = 2:3, digits = 3)
EnhancedVolcano::EnhancedVolcano(
res,
lab = res$SYMBOL,
x = "log2FoldChange",
y = "padj",
pCutoff = padj_cutoff,
FCcutoff = lfc_cutoff,
)

state_d12_vs_MEF
res <- lfc_dfs[[x]]
padj_cutoff <- 0.01
lfc_cutoff <- 0.5
res %>%
select(SYMBOL, log2FoldChange, padj) %>%
filter(padj <= padj_cutoff & abs(log2FoldChange) >= lfc_cutoff) %>%
arrange(padj) %>%
DT::datatable() %>%
DT::formatSignif(columns = 2:3, digits = 3)
EnhancedVolcano::EnhancedVolcano(
res,
lab = res$SYMBOL,
x = "log2FoldChange",
y = "padj",
pCutoff = padj_cutoff,
FCcutoff = lfc_cutoff,
)

state_iPSC_vs_MEF
res <- lfc_dfs[[x]]
padj_cutoff <- 0.01
lfc_cutoff <- 0.5
res %>%
select(SYMBOL, log2FoldChange, padj) %>%
filter(padj <= padj_cutoff & abs(log2FoldChange) >= lfc_cutoff) %>%
arrange(padj) %>%
DT::datatable() %>%
DT::formatSignif(columns = 2:3, digits = 3)
EnhancedVolcano::EnhancedVolcano(
res,
lab = res$SYMBOL,
x = "log2FoldChange",
y = "padj",
pCutoff = padj_cutoff,
FCcutoff = lfc_cutoff,
)

state_ESC_vs_MEF
res <- lfc_dfs[[x]]
padj_cutoff <- 0.01
lfc_cutoff <- 0.5
res %>%
select(SYMBOL, log2FoldChange, padj) %>%
filter(padj <= padj_cutoff & abs(log2FoldChange) >= lfc_cutoff) %>%
arrange(padj) %>%
DT::datatable() %>%
DT::formatSignif(columns = 2:3, digits = 3)
EnhancedVolcano::EnhancedVolcano(
res,
lab = res$SYMBOL,
x = "log2FoldChange",
y = "padj",
pCutoff = padj_cutoff,
FCcutoff = lfc_cutoff,
)
